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ABSTRACT: The interaction of the genes involved in intestinal development is the molecular basis of the regulatory mechanisms of 

intestinal development. The objective of this study was to identify the significant pathways and key genes that regulate intestinal 
development in Landrace piglets, and elucidate their rules of operation. The differential expression of genes related to intestinal 
development during suckling time was investigated using a porcine genome array. Time sequence profiles were analyzed for the 
differentially expressed genes to obtain significant expression profiles. Subsequently, the most significant profiles were assayed using 
Gene Ontology categories, pathway analysis, network analysis, and analysis of gene co-expression to unveil the main biological 
processes, the significant pathways, and the effective genes, respectively. In addition, quantitative real-time PGR was carried out to 
verify the reliability of the results of the analysis of the array. The results showed that more than 8000 differential expression transcripts 
were identified using microarray technology. Among the 30 significant obtained model profiles, profiles 66 and 13 were the most 
significant. Analysis of profiles 66 and 13 indicated that they were mainly involved in immunity, metabolism, and cell division or 
proliferation. Among the most effective genes in these two profiles, GN161469, which is similar to methylcrotonoyl-Coenzyme A 
carboxylase 2 (beta), and U89949.1, which encodes a folate binding protein, had a crucial influence on the co-expression network. (Key 
Words: Intestinal Development, Microarray, Gene, Expression Profile, Network) 



INTRODUCTION 

The pig (Sus scrofa) has been widely used as a model 
for studying the development of human organs and disease 
progression, owing to the similarity in size and physiology 
of the two species (Butler et al., 2007; Lunney, 2007; Ojeda 
et al., 2008). The Landrace is one of the pig breeds that are 
used world-wide in research and production because of its 
genetic stability and high lean percentage. 

The alimentary tract is essential for absorption of 
nutrients and the function of the immune system. A notable 
feature of the intestine is the complex changes that occur 
over the lifetime of an individual (Weaver et al., 1991; 
Thompson et al., 2008). The development of the intestine is 
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regulated by a large number of factors, including nutrients, 
the micro flora, the epithelium, innate and adaptive 
immunity (Pacha, 2000; Caicedo et al., 2005; Donovan, 
2006; Commare et al., 2007; Ojeda et al., 2008), which 
constitute the intestinal ecosystem. It is known that the 
adaptive change involved in enteric development is highly 
conserved and is a complex process (Gilbert and Lloyd, 
2000). The underlying molecular mechanism of intestinal 
development is associated with genes and the corresponding 
proteins that participate in this change, and with the 
interactions and signal transduction processes of the gene 
and/or protein networks. It is almost impossible. However, 
to imderstand the complex interactions of genes and 
proteins by studying them in an individual and static 
manner. The study of gene regulatory networks makes it 
feasible to analyze the interactions of genes in an integral 
and dynamic way (Tomita et al., 1999; Gilbert and Lloyd, 
2000). Gene regulatory networks have attracted 
considerable research interest owing to the rapid 
accumulation of genomic information (Kitano, 202). 
Comprehensive assessment of the gene expression profiles 
needs to be carried out in order to identify gene networks 
(Kitano, 2002), and the high-throughput techniques of 
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molecular biology can be used to provide insight into such 
networks. 

Microarray technologies have enabled researchers to 
explore the dynamics of transcription on a genome-wide 
scale (Young, 2000; Jiang and Deyholos, 2006; Zhang et al., 
2008), and they provide a systematic insight into the crucial 
factors that control development. They can be used to 
elucidate the mechanisms of intestinal development at a 
molecular level. Studies of gene expression profiling using 
microarray technologies in the development of plants and 
the human fetal pancreas have made significant progress 
(Sarkar et al., 2008; Wechter et al., 2008), but there have 
been no reports of such studies on the development of the 
porcine intestine. In the present study, we aimed using 
microarray technology to identify the significant pathways 
and key genes that regulate the intestinal development of 
Landrace piglets. 

MATERIALS AND METHODS 

Animals and tissue preparation 

All animals used in this study were managed humanely 
according to the established guidelines of the Animal Care 
Committee, Sichuan Agricultural University, P. R. China. 
The four Landrace sows were mated with a boar of the 
same breed. Two male piglets were sacrificed at 0, 3, 8, 14 
and 21 d after birth, respectively. The mid-section of the 
jejunum of each piglet was dissected rapidly. The samples 
were snap-frozen in liquid nitrogen and stored at -80°C 
until further use. 

Microarray analysis 

Isolation of RNA was performed with the QIAGEN 
RNeasy mini kit (QIAGEN, Valencia, CA, USA) following 
the manufacturer's instructions. The RNA samples were 
then analyzed using the Affymetrix GeneChip® Porcine 
Genome Array (Affymetrix, Santa Clara, CA, USA), which 
contains 23,937 probes sets that interrogate approximately 
23,256 transcripts from 20,201 genes of S. scrofa. Statistics 
were computed using the GeneChip Operating Software 
(GCOS). The data discussed in this article have been 
deposited in the NCBIs Gene Expression Omnibus (GEO, 
http://www.ncbi.nlm.nih.gov/geo/) and are accessible 
through GEO Series accession number GSE13456. 

Time sequence profile analysis of gene expression 

The differentially expressed genes from four time points 
(d 3, 8, 14, and 21) were compared with those of d 0 using 
Affymetrix GCOS software. The raw values for gene 
expression were converted to log2ratio. According to the 
analysis strategy for an experiment representing a cluster of 
short time series (Emst et al., 2005), 80 profiles were 
defined in the experiment. The profiles were analyzed that 



had the most significant probability assigned in comparison 
with that expected by the Fisher exact test and the multiple 
comparison test. If the p-value of the profiles was less than 
0.0 1, the gene expression profiles were identified to be 
significantly different with respect to time. 

Functional analysis of significant profiles of co- 
expressed genes 

According to previous studies (Yu et al., 2003; Pan et al., 

2008; Rawat et al., 2008), co-expressed genes have similar 
functions. The important functions were analyzed for those 
profiles found to be significant using Gene Ontology (GO), 
which organizes genes into hierarchical categories based on 
biological processes, molecular functions, and cellular 
components (Ashbumer et al., 2000). To undertake GO 
analysis, the p-values from the Fisher's exact test for over- 
representation of the selected genes in all GO biological 
categories are computed, as described previously. 
Significance of the GO was accepted when the p-value was 
less than 0.05. Functional annotation can provide a detailed 
description of the roles of genes located in deeper 
categories of GO, according to the hierarchical organization 
of GO. Therefore, the enrichment value was defined to filter 
the deeper categories of GO: the higher the enrichment 
value of GO, the more target GO were expected. Within a 
given category, the enrichment Re is given by equation, 

(^] 

where \^ ) nf is the number of differentially 
expressed genes within the particular category, n is the total 
number of genes within the same category, Nf is the number 
of different genes on the entire microarray, and N is the 
total number of genes on the microarray. 

Pathway analysis of differentially expressed genes 

A pathway is defined as a functional unit in which genes 
interact with each other and achieve a process of signal 
transduction. To explain test phenotypic features with 
respect to development processes, significant pathways are 
identified from differentially expressed genes. This is based 
on the p-values of the Fisher's exact test for over- 
representation of the selected genes in all pathway 
categories, as described perviously (Shalgi et al., 2007). A 
significant pathway was specified when the p-value was 
less than 0.05. To assess the significance of a particular 
category with respect to random chance, the false discovery 
rate (FDR) was estimated for a set of categories. After 
resamphng 5,000 times, the FDR was defined by equation, 
where Nk refers to the number of Fisher's test p-values less 

FDR = l--^ 
than x2 test p-values ^ . 

The p-values obtained using Fisher's exact test are 
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lower than x2 test p- values. Pathway enrichment analysis 
was performed according to a method reported previously 
(Kanehisa and Goto, 2000) to investigate whether particular 
functions were enriched among the genes that control 
distinctive characters in groups of genes that were 
differentially expressed on d 3, 8, 14, and 21, in comparison 
with d 0 (the control time point). With increasing 
enrichment, the corresponding function becomes more 
significant. Within a given category, the enrichment Re is 
given by equation, where nf is the number of different 
genes within the particular pathway category, and n is the 
total number of genes within the same category, Nf is the 
number of different genes on the entire microarray, and N is 
the total number of genes on the microarray. 

Construction and topological analysis of the gene co- 
expression network 

We transformed the normalized expression values of 
profile 66 and profile 13 from the individual correlation 
analysis into measures of pair-wise connection strengths for 
the construction of a co-expression network (Carlson et al., 
2006). The network edges were specified to be more than 
0.9, which is the value of the coefficient correlation that 
incidates a strong relationship, such as co-expression of 
genes. The centrality of a network is represented by the 
central degree (Barabasi and Oltvai, 2004). It is possible to 
determine the characteristic variables that are related to the 
distances among genes. The core of maximum order is 
referred to as the main core, or the highest k-core of the 
graph (Huber et al., 2007). A k-core subgraph of a graph can 
be generated by deleting the vertices recursively from the 
graph whose degree is less than k (Huber et al., 2007). 
Larger values of "coreness" clearly correspond to vertices 
of larger degree and a more central position in the structure 
of the network. In this work, we applied the notion of the k- 
core subgraph to predict similarity of gene function (Altaf- 
Ul-Amin et al., 2006). The nodes labeled by the same color 
may represent similar gene ontology terms. The highlight 
ontology of nodes of the same color was assessed by 
counting the genes of the same GO that had the same color. 

Real-time quantitative PCR 

Fourteen differentially expressed genes (AREG 
(Amphiregulin), ATRN (Attractin), CD59 (CD59 molecule), 
CLU (Clusterin), CST3 (Cystatin C), DPEPl (Dipeptidase), 
EPHXl (Epoxide hydrolase), LEAP2 (Liver expressed 
antimicrobial peptide 2), LOC396871 (Arginine rich 
antibacterial peptides), MS4A2 (Membrane-spanning 4- 
domains), PIAP (Putative inhibitor of apoptosis), RETN 
(Resistin), UBP (Ubiquitin-specific protease), FBP (Folate 
binding protein) dentified in the array experiment were 
selected and analyzed by RT-qPCR. Aliquots of the same 
RNA samples used for the microarray analysis were used 



for the real-time PCR. The reverse transcriptions (RT) were 
carried out with the total cellar RNA using TRIZOL 
Reagent® (Invitrogen California, USA). Each sample of 
cDNA was diluted 1:4 fold in sterile ddH20, and 1 (j1 of 
this dilution was used as the template for qPCR. Primers for 
the PCR reactions were designed by DNAStar (DNASTAR 
Inc) to have a Tm of 62 to 65 °C and an optimal annealing 
temperature of 63 °C, with a length of 100 and 250 bp. Real- 
time PCR was performed with SYBR® Premix Ex TaqTM 
(TaKaRa, Dalian, China) in 20 )il reaction volumes using 
the iCycler System (BioRad, Richmond, CA, USA) 
according to the manufacturer's instructions. Each PCR 
reaction contained 1 |li1 of cDNA, 0.2 jxM of each of the 
primers (Table 1) and 10 |il of PCR master mix. The initial 
denaturing time was 95 °C per minute, followed by 40 
cycles consisting of 94°C for 15 s, 63°C for 20 s, and 72°C 
for 30 s with a single measurement of fluorescence. 

RESULTS 

Microarray analysis 

To characterize the genes involved in intestinal 
development during the first 21 d in piglets' life, total RNA 
was extracted from the mid-section of the jejimum on post- 
natal d 0, 3, 8, 14 and 21. Gene expression was analyzed 
using a porcine genome array and statistics were computed 
using GCOS. In total, 8,427 significant transcripts were 
identified. Among these statistically significant expression 
ratios, probes representing 4,789 unique genes increased in 
transcript abundance at least 2.0 fold at one or more time 
points, while probes representing 3,618 unique genes 
decreased in abundance at least 2.0 fold. Thus, at least 24% 
of the total number of transcripts from S. scrofa were 
strongly (i.e. >2.0 fold) induced by intestinal development; 
however, at least 18% of transcripts were repressed by an 
equivalent amount (Table 2). 

Gene expression profiles significant to the development 
process 

During the course of intestinal development from d 0 to 
21, 80 model profiles were identified in which we defined 
the value of maximum unit change to be 1. For these model 
profiles, 30 were identified as significant (Figure 1). Among 
them, the two most significant profiles were profile 66 and 
profile 13, according to ascending order of the p-values 
(Figure 2). If the p-value of a profile was less than 0.01, the 
gene profile was identified as one that changed significantly 
with time. According to this principle, the lowest p-value 
profiles may represent the regulatory genes that affect the 
development of the porcine 539 genes, the expression of 
which increased significantly at d 3 and maintained a 
similar level of expression until d 21. In contrast, profile 13 
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Table 1. Primer sequences and PCR product sizes of genes selected for validation by RT 



Gene 



GenBank ID Primer sequence 



Annealing tm Product size 
(°C) (bp) 



P-actin 

Amphiregulin AREG 
Clusterin CLU 
Dipeptidase DPEPl 

Epoxide hydrolase 

EPHXl 
Liver expressed 

antimicrobial peptide 2 

LEAP2 
Arginine rich antibacterial 

peptides LOC396871 
Membrane-spanning 

4-domains MS4A2 
Putative inhibitor 

of apoptosis PIAP 
Resistin RETN 



Ubiquitin-specific 
protease UBP 
CD59 molecule CD59 

Attractin ATRN 

Folate binding protein 
FBP 
Cystatin C CST3 



SSU07786 

AY028311.1 

NM_2 1397 1.1 

NM_214108 

NM_214355 

NM_213788 

NM_214355.1 
NM_214092.1 
NM_214181 
NM_213783.1 

NM_213826.1 
NM_214170 
AB271953 
NM_213853 
NM_00 1044602 



Forwards ' -GCTGGCCGGG ACCTGACCGACT-3 ' 
Reverses ' -TTGGC ATAATTGTTAC ACGTTTGG-3 ' 
ForwardS'-GGCGCCCGTGGTGCTGTCACTC-3' 
Reverses ' -CGCTAGC AGGGGG AGCCTC ACTT-3 ' 
Forwards ' - AAGTCCCGCTTCGCCCGGAAC ATC-3 ' 
Reverses ' -CCTC AGGC ATCCCGTGGAGTTGTG-3 ' 
Forwards '-CAAC ACGCCCTGGGCCGAC AACTG-3 ' 
Reverses ' - ACGGGGGCCTGGGAC AGCTTC AGA-3 ' 
ForwardS'-CCCGTCCCCCTGGCCTACAT-3' 
Reverses '-AAAGCCGGTGGGCACGTGGACCTT-3 ' 
Forwards '-ACGGCTCCCCGATACCGGAACGAG-3 ' 
ReverseS'-TCCTGGGCCACACTTAGGGAACAGC-3' 

Forwards ' -TGCCGGAAC ATCTATGAGAATGACT-3 
Reverse 5 '-AATCAGGGGC AGAAGCTATCTCCTT-3 ' 
Forwards '-TGCCGGAACATCTATGAGAATGACT-3 ' 
Reverses ' - AATC AGGGGC AG A AGCTATCTCCTT-3 ' 
Forwards '-GCAGCCCGCTTTAAAACATTCTGTA-3 ' 
Reverse S '-TTGAACGCGACTGATGAACTCCTGT-3 ' 
Forwards ' -GGGGCTGCTGGTGTGGGGC AAGT-3 ' 
Reverses ' -CGCGC ACGTCCC AAGAGCC AC AG-3 ' 
ForwardS'-GGGCCCGGTCGGTTTGCACAACAT-3' 
Reverses ' -CACCGCCTTCTGCCGGCTGTCCT-3 ' 
Forwards ' - AAGCCGTGCC ACCC AAAACTTACTA-3 ' 
Reverses ' -C AGC AGAGCGGTTTTCCCTGATGA-3 ' 
Forward S ' -TGGCC AGTCGTCCGTTTGCCTCTG-3 ' 
Reverse S ' -AGGGCGCTGGCC ACTGCGAGAC-3 ' 
Forward: S'-ACGGGGACCACCTGCGACACATT-3' 
Reverse: S ' -CGC ACTC AAGGCC AAGGCGTAGA-3 ' 
Forward: S '-CCCAACAGGCAAGGCAGCTAAAGC-3 ' 
Reverse: S ' -CTGCC AAGACCGCGTGGTTGAGATT-3 ' 



63.5 
61.5 
61.4 
64.1 
62.3 
59.9 

57 
54.8 
56.1 
63 

61 

56.2 

61.5 

66 

59.3 



154 
179 
225 
183 
243 
169 

250 
244 
240 
218 

199 
167 
220 
1,320 

203 



had a significance level of 8.8x10" and contained 351 
genes, the expression of which decreased significantly at d 
3 and maintained a similar level of expression until d 21. 

Significant gene functions and genes of significant 
profiles 

Results of the analysis of GO categories for profile 66 
implied that the main enriched GO categories were focused 



on control checkpoints for cell size, cell wall biogenesis, 
deoxyribonucleotide catabolic processes, fungal-type cell 
wall biogenesis, fungal-type cell wall organization and 
biogenesis, lymphocyte chemotaxis, M phase-specific 
microtubule processes, membrane depolarization, 
microvillus biogenesis, mitochondrial depolarization, 
mitotic metaphase, etc. during growth ftom 0 to 21 d (Table 
3). To further investigate the most significant functions, the 



Table 2. Summary of data obtained from the ] 


genechip analysis 










d 3 vs d 0 


d 8 vs d 0 


d 14 vs d 0 


d21 vsdO 


Total transcripts 


20,201 


20,201 


20,201 


20,201 


Differentially expressed genes (DEG) 


1,881 


1,918 


2,204 


2,404 


DEG/total (%) 


9.31 


9.49 


10.91 


11.9 


Increased genes (IG) 


1,117 


994 


1340 


1338 


IG/DEG (%) 


59.38 


51.82 


60.8 


55.66 


Decreased genes (DG) 


764 


924 


864 


1,066 


DG/DEG (%) 


40.62 


48.18 


39.2 


44.34 
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Profiles ordered based on the p-value significance of number of genes assigned versus expected 
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Figure 1. Model expression profiles of differentially expressed genes. The number in the top left hand comer of a profile box is the 
model profile ID number, which has no special meaning; the lower left hand comer contains the p-value. Model expression profiles are 
ordered from left to right and from top to bottom according to gradually increasing significance. Significant profiles have statistical 
significance (corrected with the multiple permutation test) less than 0.01. 



annotation of the GO was analyzed using Amigo (a website 
tool (http ://aniigo . geneontology.org/cgi-bin/amigo/go .cgi) . 
It was concluded that the process of the cell cycle initiated 
at a checkpoint and the immune response to fungi 
stimulated by fungal infection were the most significant 
phenomena that occurred from d 3 to d 2 1 . 

In contrast, the main enriched GO categories for profile 
13 were focused on folic acid transport, negative regulation 
of endocytosis, the phosphatidylcholine biosynthetic 
process, phosphatidylcholine metabolic processes, the 
response to cobalt, manganese and zinc ions, catabolic 
processes for aminoglycans, glycosaminoglycans, and 
leucine, leucine metabolic processes, a lipopolysaccharide- 
mediated signaling pathway, Mo-molybdopterin cofactor 
biosynthetic processes, and Mo-molybdopterin cofactor 
metabolic processes, etc., as presented in Table 4. 

Distribution of significant pathways of differentially 

expressed genes 

The pathway was analyzed using the KEGG and 
Genmapp databases. The significant pathways for 
differentially expressed genes were applied to the Fisher 
exact test and tests of multiple comparisons to obtain the 
significant functions. From the results of Table 6, it can be 



seen that the metabolism of xenobiotics by cytochrome 
P450, Alzheimer's disease, antigen processing and 
presentation, the PPAR signaling pathway, complement and 
coagulation cascades, cell adhesion molecules (CAMs), 
autoimmune thyroid disease, allograft rejection, graft- 
versus-host disease, natural killer cell mediated cytotoxicity, 
and cytokine-cytokine receptor interaction were the 
significant pathway categories. It is well known that owing 
to the changes in diet and living environment that occur 
after the birth of piglets, adaptive changes occur in the 
internal environment. In this process various immune 
reactions are activated, such as antigen processing, 
complement cascades, allograft rejection, graft- versus-host 
disease, natural killer cell mediated cytotoxicity and 
cytokine-cytokine receptor interaction. It can be concluded. 
Therefore, that the immune response may be coupled with 
the growth of the organism (Table 5). 

The gene co-expression network 

The gene co-expression network was constructed with 
respect to gene function associations. Given that the 
elements of a network represent a variety of gene regulatory 
abiUties, a large-scale gene network can be divided into 
subgraphs, named k-core networks, in which all genes 



1486 



Sun etal. (2012) Asian-Aust. J. Anim. Sci. 25:1481-1492 




Figure 2. Gene clusters of significant expression profiles 66 and 13. The horizontal axis represents time, and the vertical axis shows the 
time series of the level of gene expression after log normal transformation. Each diagram contains the model profile ID at the top center 
"Genes assigned" represents the number of genes actually assigned to the profile; "Genes expected" represents the expected number of 
genes based on a permutation test of time points; p-value represents the significance of the number of genes assigned versus expected. 



affect each other (Huber et al., 2007). Generally, the 
complexity of gene relationships increases with the rank of 
the k-core value. In the light of the definition of the k-core 
network, the core status within a large-scale gene network 
consists of subgraphs of higher k value. In this work, we 
used the k-cores of protein-protein interaction networks to 



define the main gene functions in the main subgraph. 
Moreover, we aimed to find the main GO assigned by the 
maximum number of genes in a separate k-core and then to 
define the key gene functions at each level of the network. 

It has been reported that the core functions at the core 
status of a network have a top k-core level (Altaf-Ul-Amin 
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Table 3. Significant GO of profile 66 



Go ontology 



~ • 

p- value Enrichment Go ontology 



Table 4. Significant GO of profile 13 



p-value Enrichment 



Cell size control checkpoint 0 5.94 

Cell wall biogenesis 0 5.94 

Deoxyribonucleotide cataboUc process 0 5.94 

Fungal-type cell wall biogenesis 0 5.94 

Fungal-tj'pe cell wall organization 0 5.94 
and biogenesis 

Lymphocyte chemotaxis 0 5.94 

M phase specific microtubule process 0 5.94 

Membrane depolarization 0 5.94 

Microvillus biogenesis 0 5.94 

Mitochondrial depolarization 0 5.94 

Mitotic metaphase 0 5.94 

Negative regulation of 0 5.94 

calcium-dependent cell-cell adhesion 

Negative regulation of 0 5.94 

megakaryocyte differentiation 

Negative regulation of 0 5.94 

mitochondrial depolarization 

Phosphatidic acid biosynthetic process 0 5.94 

Phosphatidic acid metabolic process 0 5.94 

Pyrimidine deoxyribonucleotide 0 5.94 

catabolic process 

Regulation of calcium-dependent 0 5.94 

cell-cell adhesion 

Regulation of 0 5.94 

mitochondrial depolarization 

Regulation of mitochondrial 0 5.94 

membrane potential 

Release of cytochrome c 0 5.94 

from mitochondria 

Deoxyribonucleoside monophosphate 0 3.96 

biosynthetic process 

Deox5Tibonucleoside monophosphate 0 3.96 

metabolic process 

Leukocyte mediated cytotoxicity 0 3.96 

Negative regulation of leukocyte 0 3.96 

mediated cytotoxicity 

Ncgativcrcgulation of T cell 0 3.96 

mediated cytotoxicity 

Regulation of cell kilUng 0 3.96 

et al., 2006), and a higher degree but lower clustering 
coefficient (Barabasi and Oltvai, 2004). A k-core 
subnetwork with a k- value of 12 was located at the core 
status within the large scale gene network (Figure 3). Using 
the hierarchical categories of gene ontology, it could be 
concluded that most of the genes in this subnetwork are 
attributed to cell communication, which may be related to 
cell division or proliferation through cell adhesion and 
transport of nutrients (for details of the genes. Tables 3 and 
4). In addition, most genes in the secondary network (with a 
k-core value of 11) play an important role in metabolic 



Phosphatidylcholine biosynthetic process 
Folic acid transport 
Response to manganese ion 
Response to cobalt ion 
Negative regulation of endocytosis 
Response to zinc ion 
PhosphatidylchoUne metabolic process 
Negative regulation of 

chemokine biosynthetic process 
Aminoglycan catabolic process 
Lipopoly saccharide- mediated 

signaling pathway 
Glycosaminoglycan catabolic process 
Regulation of lipopolysaccharide- 

mediated signaling pathway 
Vitamin transport 
Molybdopterin cefaclor 

biosynthetic process 
Leucine catabolic process 
Mo-molybdopterin cofactor 

metabolic process 
Negative regulation of lipopolysaccharide- 

mediated signaling pathway 
Negative regulation of INK activity 
Molybdopterin cofactor metabolic process 
Leucine metabolic process 
Mo-molybdopterin cofactor 
biosynthetic process 
Cellular zinc ion homeostasis 



5.94 
5.94 
5.94 
5.94 
5.94 
5.94 
5.94 
2.97 

2.97 
2.97 

2.97 
2.97 

2.97 
2.97 

2.97 
2.97 

2.97 

2.97 
2.97 
2.97 
2.97 

2.38 



processes of nitrogen compounds and amines (for detailed 
information. Tables 3 and 4), except for those genes without 
detailed functional annotation. Among these genes, 
CN161469 is similar to methylcrotonoyl-Coenzyme A 
carboxylase 2 (beta) when analyzed by BLAST, and it takes 
part in the regulation of other genes within the secondary 
level subnetwork because of the higher value of its degree. 
Most notably, the gene U89949.1, which encodes a folate 
binding protein, possessed the highest k-core and degree, 
and the lowest clustering coefficient (Table 6). These 
features indicate that this is the core gene that is located at 
the center of both the large scale network and the 12 k-core 
subnetwork, and it regulates directly 28 neighboring genes. 

VaUdation of the gene chip data by RT-PCR 

To confirm whether the genes identified were 
differentially expressed, quantitative RT-PCR analysis for 
specific transcripts was carried out as described previously. 
Based on our microarray analysis, 14 genes, including three 
antimicrobial peptides (LEAP2, LOC396871, and MS4A2), 
two proteins involved in the regulation of proliferation and 
apoptosis of epithelial cells (PIAP and AREG), and genes 
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Table 5. Significant pathway categories (p<0.01, FDR<0.05) 



Pathway p-value Enrichment FDR 



MetaboUsm of xenobiotics 


0 


4.25 


0 


by cytochrome P450 








Alzheimer's disease 


0 


4.12 


0 


Antigen processing 


0 


4.1 


0 


and presentation 








PPAR signaling pathway 


0 


3.89 


0 


Complement and 


0 


3.57 


0 


coagulation cascades 








Cell adhesion molecules 


0 


3.47 


0 


(CAMs) 








Autoimmune thyroid disease 


0 


3.36 


0 


Allograft rejection 


0 


3.36 


0 


Graft-versus-host disease 


0 


3.27 


0 


Natural killer cell 


0 


2.7 


0 


mediated cytotoxicity 








Cytokine-cytokine 


0 


2.19 


0.01 


receptor interaction 









associated with a number of physiological and metabolic 
processes (ATRN, CD59, CLU, DPEPl, EPHXl, RETN 



and UBP) were selected. The core gene FBP was also 
involved in the RT-PCR analysis. For all of these genes, the 
expression ratios measured by RT-qPCR and by microarray 
were highly correlated (r = 0.74, 0.73, 0.95, 0.93 at d 3, 8, 
14, 21 respectively) (Table 7). Genes involved in apoptosis, 
signal transduction, and energy and protein metabolism, 
such as CLU, MS4A2 and UBP, were up-regulated during 
suckling, as indicated by both the Genechip and RT-qPCR 
experiments. However, HCRTR2 displayed significant 
differential expression when detected by RT-qPCR but no 
significance on Genechip testing. The correlation of the 
results from the RT-qPCR and Genechip testing helped us to 
analyze the reliability of the assays so as to reveal true 
differences in the gene expression profiles in the intestine. 

DISCUSSION 

Owing to the continuous improvement of microarray 
technology, the rapid development of bioinformatics and the 
depth of research on genomics and proteomics, global 
genome expression profile analysis provides a more 



Table 6. Nclwork slruclurc paraniclcr lists for each gene 



Probe Set ID 


NCBl-gcnelD 


Clustering coefficient 


Degree 


k-core 


Ssc. 14544. l.Sl_at 


U89949.1 


0.58 


28 


12 


Ssc.21579.1.Sl_at 


CF789025 


0.58 


27 


12 


Ssc.l735.1.Sl_at 


Bill 8904 


0.62 


26 


12 


Ssc.5250.1.Sl_at 


BQ604261 


0.58 


26 


12 


Ssc.l5069.1.Sl_at 


BI345060 


0.63 


25 


12 


Ssc.22910.1.Al_at 


C0986224 


0.67 


24 


12 


Ssc. 10748. l.Al_at 


BQ597862 


0.61 


24 


12 


Ssc.l6817.1.Sl_at 


BE233328 


0.64 


23 


12 


Ssc.25002.1.Sl_at 


BG833819 


0.69 


22 


12 


Ssc.6080.1.Sl_at 


NM_213870.1 


0.68 


22 


12 


Ssc.l5430.1.Al_at 


CO937301 


0.67 


22 


12 


Ssc. 1 1246. l.Al_at 


BI181165 


0.55 


22 


12 


Ssc.l8215.1.Sl_at 


NM_214365.1 


0.71 


21 


12 


Ssc.l3713.1.Al_at 


BQ603271 


0.62 


20 


12 


Ssc.3714.1.Sl_a_at 


CK456589 


0.71 


19 


12 


Ssc.21141.1.Sl_at 


BX666500 


0.6 


19 


12 


Ssc.26028.1.Sl_at 


BX926563 


0.43 


19 


9 


Ssc.l209.1.Sl_at 


CK465976 


0.76 


18 


12 


Ssc.l8033.1.Al_at 


CFl 80565 


0.75 


18 


12 


Ssc.5052.1.Sl_at 


BI398918 


0.69 


18 


12 


Ssc.2795.1.Sl_at 


BI183754 


0.68 


18 


12 


Ssc.21952.1.Sl_at 


BX666921 


0.54 


18 


10 


Ssc.9503.1.Al_at 


CF789586 


0.46 


18 


9 


Ssc.24263.2.Sl_at 


AW480085 


0.55 


16 


9 


Ssc.l2622.2.Al_a_at 


CO948012 


0.52 


16 


9 


Ssc. 27572. l.Sl_at 


CN161469 


0.67 


15 


11 


Ssc. 1 0988. 1.Sl_al 


C0947173 


0.72 


15 


10 



Representative public ID can be queried on NCBI. Clustering coefficient and degree describe the network properties of genes. K-core represents the 
subnetwork in which genes are located. 
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Figure 3. Gene association network. Nodes denote genes; undirected links represent gene-gene inten'elation; the size of the nodes 
represents the power of the interrelationship among the nodes. Different colors of nodes represent distinct k-core gene classes. 



Table 7. Genes differentially expressed in the Genechip data 



Symbol 


Method 


0-d 


3-d old 


8-d old 


14-d old 


21-d old 


Clusterin 


RT-PCR 




0.584 


0.822 


8.830 


10.255 


CLU 


Genechip 




1.546 


2.923 


8.558 


13.860 


Amphiregulin 


RT-PCR 




2.263 


1.004 


1.546 


1.969 


AREG 


Genechip 




10.317 


1.950 


5.347 


7.130 


Folate binding 


RT-PCR 




1.005 


0.380 


0.850 


3.235 


FBP 


Genechip 




0.295 


0.304 


0.264 


0.2606 


Dipeptidase 


RT-PCR 




0.115 


0.111 


0.179 


0.251 


DPEPl 


Genechip 




0.255 


0.248 


0.375 


0.468 


Membrane-spanning 4-domains 


RT-PCR 




0.475 


2.670 


31.138 


7.704 


MS4A2 


Genechip 




0.831 


1.798 


8.877 


4.304 


Utative inhibitor of apoptosis 


RT-PCR 




0.905 


1.449 


6.358 


2.027 


PIAP 


Genechip 




2.366 


2.106 


4.854 


3.761 


Ubiquitin-specific protease 


RT-PCR 




5.199 


1.770 


5.554 


28.016 


UBP 


Genechip 




4.774 


3.007 


3.499 


16.768 


Attractin 


RT-PCR 




0.211 


0.477 


0.758 


0.197 


ATRN 


Genechip 




0.239 


0.679 


0.186 


0.361 


Liver expressed antimicrobial peptide 2 


RT-PCT 




1.775 


4.712 


2.652 


2.841 


LEAP2 


Genechip 




3.882 


3.552 


2.484 


2.758 


Resistin 


RT-PCT 




0.995 


2.634 


1.1760 


0.309 


RETN 


Genechip 




2.874 


3.181 


3.077 


2.113 


CD59 molecule 


RT-PCT 




0.080 


0.099 


0.346 


0.197 


CD59 


Genechip 




1.485 


0.848 


1.048 


0.992 


Epoxide hydrolase EPHXl 


RT-PCT 




6.169 


22.177 


64.977 


33.188 




Genechip 




12.195 


26.257 


25.856 


35.314 


Ai'ginine rich antibacterial peptides 


RT-PCT 




4.342 


23.007 


10.273 


1.649 


LOC396871 


Genechip 




3.211 


3.927 


4.822 


1.118 


Cystatin C 


RT-PCR 




1.588 


0.718 


2.674 


8.490 


CST3 


Genechip 




1.474 


2.631 


2.469 


1.996 


r 






0.74 


0.72 


0.95 


0.93 
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comprehensive and dynamic analysis of intestinal 
development at the genome-wide transcriptional level 
(Ravasz et al., 2002; Schlitt et al., 2003; Stears et al., 2003; 
Barabasi and Oltvai, 2004; Han, 2008; Lee et al., 2008). 
The understanding and description of cells, organs and 
living individuals on a system level is a new area of 
biological research (Davidson et al., 2002; Kitano, 2002; 
Oltvai and Barabasi, 2002). Various kinds of model 
organism have been proved to be useful in this area 
(Covington et al., 2008; Sarkar et al., 2008; Schweikl et al., 
2008; Wechter et al., 2008). To our knowledge, this is the 
first study involving whole genome analysis of the genes 
related to intestinal development in Landrace piglets at 
different stages of growth (0, 3, 8, 14, and 21 d after birth). 

In the present study, more than 8,000 differentially 
expressed transcripts were identified by microarray 
technology. This technique was validated by RT-qPCR for 
utilization in a comprehensive study of gene expression 
profiles. However, many of these transcripts had no detailed 
information in the cDNA resources deposited for swine, 
which limited the analysis of their biological importance. 
This indicates that many porcine genes have yet to be 
analyzed. Complex molecular events pertinent to intestinal 
development occur during the suckling stage in piglets 
(Hall and Byrne, 1989; Buddington, 1994; Thaler and 
Cummings, 2008). We analyzed the development process of 
the gene expression profile, and 80 model profiles were 
obtained in all. These profiles do not cover all the possible 
profiles that could be generated in the intestine of Landrace 
piglets during suckling. Among these profiles, 30 were 
identified as significant and two of them (profiles 66 and 1 3) 
were highly significant. Significant profiles may indicate 
that common fiinctions were mainly attributed to the co- 
expressed genes, and these functions influence the 
biological characteristics of the organism (Gracey et al., 
2004). We can therefore conclude, according to the 
significance of the profiles, that the genes assigned to these 
two profiles may be highly specific for the intestinal 
development of piglets in early life. In order to obtain a 
better understanding of the genes assigned to these two 
profiles and to elucidate the molecular mechanism of 
development of the intestine, analysis of GO categories, 
pathway analysis, and topological analysis of the gene co- 
expression network were carried out. 

Annotation of GO categories has proved to be 
remarkably useful for the mining of functional and 
biological significance from the results of microarray 
studies. In our study, significant GO categories were 
determined by the hypergeometric distribution and multiple 
comparison tests. The results obtained for the GO 
categories implied that the genes in profiles 66 and 13 were 
involved in many biological processes, among which the 



distinguished phenomena were initiation of the cell cycle at 
a checkpoint and the immune response to fungi for profile 
66, while the folic acid transport process was most notable 
for profile 13. These results suggest that the genes in profile 
66 may regulate cell proUferation and differentiation, and 
the genes in profile 13 may modulate the metabolism of 
nutrients. 

Analysis of the pathway for the lists of differentially 
expressed genes may allow us to find the target pathway 
that regulates phenotypic differences. Pathway analysis in 
the present study indicated that the genes in profiles 66 and 
13 focused mainly on the immune system, metabolism, cell 
adhesion molecules, etc. It may be concluded that these 
genes play an important role in the above pathways, which 
agrees with the results obtained from the analysis of the GO 
categories. 

During biological processes, a macromolecular network 
can be constructed according to the yeast two-hybrid 
method (Y-2H) (Smidtas et al., 2006) or an algorithmic 
prediction based on the gene function correlation and 
expression profiles (Nikiforova and Willmitzer, 2007). The 
complexity of the network model based on the high 
throughput gene expression test using the algorithmic 
prediction made it feasible to determine snapshots of 
protein-protein interactions, gene expression regulatory 
networks, and metabolic networks among different gene 
groups (Carlson et al., 2006). To find the most effective 
genes in profiles 66 and 13, the gene co-expression network 
was constructed using the algorithmic prediction, and 
topological analysis was carried out, the results of which 
demonstrated that more than 20 effective genes were 
identified. Most of these genes were attributed to cell 
communication, which may be related to cell division or 
proliferation through cell adhesion and the transport of 
nutritients. Among these genes, CN161469 was found to be 
similar to methylcrotonoyl-Coenzyme A carboxylase 2 
(beta), which suggests that it possesses the similar 
biological functions. In addition, the gene U89949.1, which 
encodes a folate binding protein that is involved in the 
transport of foliate (Vallet et al., 2001; Kim and Vallet, 
2004), played an major role in the network. This finding 
was in accordance with the results of the assessment of the 
GO categories and pathway analysis. It has been 
demonstrated that down-regulated genes reflect a state of 
malnutrition, which may be induced by a lack of the 
essential molecules for growth, such as folic acid, trace 
elements, and aminoglycan, in early life. In profile 13, the 
gene U89949.1 was down-regulated, which may imply a 
decrease in foUc acid transport that may subsequently 
induce malnutrition. This gene therefore provides a channel 
that allows regulation of the metabolism on a molecular 
level to balance malnutrition and over-nutrition. 
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IMPLICATIONS 

In order to increase our understanding of the molecular 
mechanisms involved in intestinal development, it is 
essential to obtain more information on the functions of 
genes in the mammalian genome and their corresponding 
products. In the present study, we have provided a rich new 
information resource, which describes the significant 
pathways and effective transcripts that are related intimately 
with porcine intestinal development. Although there was a 
certain degree of ambiguity when we established the 
regulatory network based on the expression profiling data, 
this represents a valuable attempt to guide further in-depth 
research into the candidate genes and signal transduction 
pathways involved in intestinal development. 
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